/* libwatopt++ - Optimization Toolkit
 *
 * Copyright (c) 2007 Christopher Alexander Watford
 * <christopher.watford@gmail.com>
 *
 * Permission is hereby granted, free of charge, to any person obtaining a copy
 * of this software and associated documentation files (the "Software"), to
 * deal in the Software without restriction, including without limitation the
 * rights to use, copy, modify, merge, publish, distribute, sublicense, and/or
 * sell copies of the Software, and to permit persons to whom the Software is
 * furnished to do so, subject to the following conditions:
 *
 * The above copyright notice and this permission notice shall be included in
 * all copies or substantial portions of the Software.
 *
 * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
 * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
 * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
 * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
 * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
 * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS
 * IN THE SOFTWARE.
 *
 * $Id: hopfield_network.cpp 2 2007-05-06 18:01:15Z christopher.watford $
 */

#include <cstdlib>
#include <ctime>
#include <cmath>
#include <limits>

#include <assert.h>

#include <hopfield_network.hpp>

#include "utils.hpp"

namespace watopt
{

   HopfieldNetwork*
   SetupKofN(size_t k, size_t n, double** II)
   {
	   assert(k <= n);
	   assert(II != NULL);

	   seed_random();

	   HopfieldNetwork* network = new HopfieldNetwork(n);
	   (*II) = network->CreateInputArray();
	   double *uu = network->Get_u(),
		   *TT = network->Get_T();
	   for(size_t ii = 0U; ii < n; ii++)
	   {
		   uu[ii] = ((2.0 * k) - n)*(network->alpha()/(double)n)
					   + ((rand01()-0.5)/100.0);
		   (*II)[ii] = ((2.0 * k) - 1.0);
		   for(size_t jj = 0U; jj < n; jj++)
		   {
			   TT[Get_T_Index(ii,jj,n)] = (ii == jj ? 0.0 : -2.0);
		   }
	   }

	   return network;
   }

   static double
   piecewise_linear(double val, double alpha)
   {
	   return (val >= alpha ? 1.0 : (val <= -alpha ? 0.0 : (val+alpha)/(2.0*alpha) ));
   }

   double*
   HopfieldNetwork::Run(double *II, size_t input_count, size_t max_iterations)
   {
	   assert(II != NULL);
	   assert(input_count == this->_nodes);

	   // initialization
	   size_t* order = new size_t[this->_nodes];
	   for(size_t aa = 0U; aa < this->_nodes; aa++)
	   {
		   order[aa] = aa;
		   _du[aa] = 0.0;
		   _VV[aa] = piecewise_linear(_uu[aa], this->_alpha);
	   }
   	
      double sumOfAbsolutedelta_u = std::numeric_limits<double>::max();
       for(unsigned int cnt = 0U; ((cnt < max_iterations) && (sumOfAbsolutedelta_u > this->_eps)); cnt++)
       {
		   // equilibrium check
 		   sumOfAbsolutedelta_u = 0.0;

		   // present the inputs in random order
		   shuffle<size_t>(order, this->_nodes);

		   for(size_t jj = 0U, ii; jj < this->_nodes; jj++)
		   {
			   ii = order[jj];

			   double feedback = 0.0;
			   for (size_t mm = 0U; mm < this->_nodes; mm++)
			     feedback += _TT[Get_T_Index(ii,mm,this->_nodes)]*_VV[mm];

			   _du[ii] = (-_uu[ii]/_tau + feedback + II[ii]) * _dt;

			   sumOfAbsolutedelta_u += std::abs(_du[ii]);

			   //Asynchrous updates
			   _uu[ii] += _du[ii];
			   _VV[ii] = piecewise_linear(_uu[ii], _alpha);
		   }
       }

	   // make a copy of V to return
	   double* vv_ret = new double[this->_nodes];
	   for(size_t kk = 0U; kk < this->_nodes; kk++)
		   vv_ret[kk] = this->_VV[kk];

	   return vv_ret;
   }

   double*
   HopfieldNetwork::CreateInputArray(void) const
   {
	   return new double[this->_nodes];
   }

   double*
   HopfieldNetwork::Get_u(void) const
   {
	   return this->_uu;
   }

   double*
   HopfieldNetwork::Get_V(void) const
   {
	   return this->_VV;
   }

   double*
   HopfieldNetwork::Get_T(void) const
   {
	   return this->_TT;
   }

   HopfieldNetwork::HopfieldNetwork(
	   size_t nodes,
	   double tau,
	   double dt,
	   double alpha,
	   double eps
	   )
   :	_nodes(nodes)
   ,	_TT(new double[nodes*nodes])
   ,	_VV(new double[nodes])
   ,	_uu(new double[nodes])
   ,	_du(new double[nodes])
   ,	_tau(tau)
   ,	_dt(dt)
   ,	_alpha(alpha)
   ,	_eps(eps)
   {
   }

   HopfieldNetwork::~HopfieldNetwork(void)
   {
	   delete _TT;
	   delete _VV;
	   delete _uu;
	   delete _du;

	   _TT = NULL;
	   _VV = NULL;
	   _uu = NULL;
	   _du = NULL;
   }

}
